clear
clc
close all
load('3.mat')
pic = double(allData1);
show_pic_demo = pic(:,:,55);
imagesc(show_pic_demo);
axis tight
axis equal
[y,x] = ginput(6);
disp(y);
xx = round(x');
yy = round(y');
[x_all, y_all, k_all] = size(allData1);

[xx_all, yy_all] = size(xx);disp(xx_all);
spectrum_mat = zeros(xx_all, yy_all, k_all);

% for ii = 1 : xx_all
%     for jj = 1 : yy_all
%         for i = 1 : x_all
%             for j = 1 : y_all
%                 if(norm([j,i] - [yy(ii,jj),xx(ii,jj)])<=4)
%                  % if(norm([j, min(3*i, 0)] - [yy(ii, jj), xx(ii, jj)]) <= 4)
%                     for kk = 1 : k_all
% 
%                         spectrum_mat(ii,jj,kk) = spectrum_mat(ii,jj,kk)+pic(i,j,kk);
%                     end
%                     show_pic_demo(i,j) = nan;
%                 end
%             end
%         end
%     end
% end
half_width = 8;
for ii = 1 : xx_all
    for jj = 1 : yy_all
        % 确定矩形的边界
        x_min = max(1, xx(ii, jj) - 1); % 上边界
        x_max = min(x_all, xx(ii, jj) + 2); % 
        y_min = max(1, yy(ii, jj) - half_width); % 左边界
        y_max = min(y_all, yy(ii, jj) + half_width); % 

        % 遍历所有像素点，累加位于矩形内的像素值
        for i = x_min : x_max
            for j = y_min : y_max
                for kk = 1 : k_all
                    spectrum_mat(ii, jj, kk) = spectrum_mat(ii, jj, kk) + pic(i, j, kk);
                end
                show_pic_demo(i, j) = nan; % 将已累加的像素设置为 NaN 以便显示
            end
        end
    end
end

imagesc(show_pic_demo);
for point_index = 1 : size(x,1)
    text(y(point_index), x(point_index), num2str(point_index), 'Color', 'red', 'FontSize', 12);
end
axis tight
axis equal
                    



index = 1;
for p = 1 : xx_all
    for q = 1 : yy_all
        s(index,:) = squeeze(spectrum_mat(p,q,:));
        index = index + 1;
    end
end

figure('Units', 'normalized', 'Position', [0.1, 0.1, 0.8, 0.4]); % 设置窗口比例
peak_range = [40, 52];
max_per_row = 5;  % 每行最多显示5个图形
total_plots = size(x, 1);  % 总共有多少个折线图要显示
% 计算总共需要多少行
num_rows = ceil(total_plots / max_per_row); 
for t = 1:size(x, 1)
    row = ceil(t / max_per_row);  % 计算当前行
    col = mod(t-1, max_per_row) + 1;  % 计算当前列
    subplot(num_rows, max_per_row, (row-1)*max_per_row + col);
    % 数据归一化处理
    s_norm(t, :) = s(t, :) ./ (max(max(s(t, 7:end))));
    
    % 找到指定范围内的峰值
   
    [nip_idx,nip] = FindPeak(s_norm(t, :), peak_range); % 找到范围内的峰值和索引
    
    % 绘图
    plot(s_norm(t, :), 'linewidth', 2);
    
    % 在指定范围内标记峰值
    hold on;
    peak_global_idx = nip_idx ; % 转换为全局索引
    plot(peak_global_idx, nip, 'ro', 'MarkerSize', 8, 'LineWidth', 2); % 标记峰值
    
    % 添加标题和标签
    title(['t = ', num2str(t), ', 峰 = ', num2str(nip_idx)]); 
    xlabel('Index'); 
    ylabel('Normalized Value'); 
    
    % 添加网格
    grid on;
    
    % 暂停，等待用户交互
    pause;
end
% for t = 1 : size(x,1)
%     subplot( 1,size(x, 1), t); % 根据 t 创建子图
%     s_norm(t, :) = s(t, :) ./ (max(max(s(t, 7:end))));
%     plot(s_norm(t, :), 'linewidth', 2);  %丢信息
%     title(['t = ', num2str(t)]); % 给每个子图添加标题
%     xlabel('Index'); % 添加横轴标签
%     ylabel('Normalized Value'); % 添加纵轴标签
%     grid on; % 添加网格
%     pause; % 暂停，等待用户交互
%     % s_norm(t,:) = s(t,:)./(max(max(s(t,7:end))));
%     % plot(s_norm(t,:),'linewidth',2);
%     % hold on;
%     % pause    
% end
